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Abstract 

We consider linear acoustic propagation at a fixed frequency in a subsonic convective flow 
around a scattering object, under the assumption that the flow is uniform far away from the 
object. Using the Prandtl-Glauert transformation, the convected Helmholtz equation yields the 
classical Helmholtz equation in the exterior domain where the flow is uniform. We then derive two 
coupled methods to solve the problem approximately using the finite clement method in the interior 
domain around the scattering object and the boundary element method on the coupling surface 
with the exterior domain. The first one admits infinitely many solutions at resonant frequencies. 
The second one is based on the combined field integral equations method, and is well-posed at 
all frequencies. We perform the analysis of the two formulations and present numerical results 
illustrating the advantages of the second formulation. 

1 Introduction 

The scope of the present work is the computation of linear acoustic wave propagation at a fixed 
frequency in the presence of a flow. The simplest model is the Helmholtz equation, when the flow is at 
rest. When the medium of propagation is not at rest, the governing equation is a convected Helmholtz 
equation resulting from the linearized harmonic Euler equations. The numerical method we propose 
hinges on the linearity of the convected Helmholtz equation. Therefore, nonlinear interaction between 
acoustics and fluid mechanics phenomena is not considered; we refer to the early work of Lighthill for 
aerodynamically generated acoustic sources [29, 30J, to [21 J for a review on nonlinear acoustics, and 
to |38| for the coupling of Computational Aero Acoustic (CAA) and Computational Fluid Dynamics 
(CFD) solvers. 

Our main assumption on the flow is that it is irrotational close to the scattering object and 
uniform far away from it. This geometric setup leads to a partition of the unbounded medium of 
propagation into two subdomains, the bounded interior domain near the scattering object where the 
flow is non-uniform and the unbounded exterior domain far away from the object where the flow 
is uniform. A convenient tool to handle the convected Helmholtz equation is the Prandtl-Glauert 
transformation, which was used in |14j in the case of globally uniform and subsonic flow. Our first 
contribution is the extension of this transformation to the case of non-uniform flow in the interior 
domain, still under the assumption that the convective flow is everywhere subsonic. In particular, 
the transformed convected Helmholtz equation simplifies into the classical Helmholtz equation in 
the exterior domain, and in the interior domain, leads to an anisotropic definite positive diffusive 
matrix and a skew-symmetric first-order term. Our second contribution is to derive and analyze a 
mathematically well-posed numerical method to approximate the transformed convected Helmholtz 
equation uniformly at all frequencies. Numerical methods for acoustic scattering include infinite finite 
elements [iQl [3] , the method of fundamental solutions [17] , and the retarded potential method [151 ES] . 



In the case of the Helmholtz equation, the problem of acoustic scattering by a solid object can be 
reduced to finding unknown functions defined on its surface and which solve integral equations. The 
boundary element method (BEM) is then used to approximate these equations numerically |35j . When 
the medium of propagation is non-uniform or when there is a non-uniform flow around the scattering 
object, a volumic resolution has to be considered, e.g. using a finite element method (FEM). If such 
non-uniformities occur only in a given bounded domain, it is possible to benefit from the advantages 
of both a volumic resolution and an integral equation formulation. Coupling BEM and FEM at 
the boundary of the given bounded domain allows this. Coupled BEM-FEM can be traced back to 
Zienkiewicz, Kelly and Bettess in 1977 [3], and Johnson and Nedelec in 1980 [25]. These methods 
have been applied to Maxwell equations |28| . In what follows, we first derive an unstable coupled 
formulation which is well-posed except at resonant frequencies for which it leads to infinitely many 
solutions. All these solutions deliver the same acoustic potential in the exterior domain, but because 
of the lack of uniqueness, the numerical procedure becomes ill-conditioned. This problem has been 
tackled in 0[22], where a stabilization of the coupling, based on combined field integral equations, has 
been proposed in the zero flow case by means of introducing an additional unknown at the coupling 
surface. Using a similar stabilization method, we derive a stable coupled formulation for the convected 
Helmholtz equation. This formulation is well-posed at all frequencies. Our numerical results show 
that the unstable formulation yields results polluted by spurious oscillations in the close vicinity of 
resonant frequencies, whereas the stable formulation yields consistent solutions at all frequencies. This 
advantage of the stable formulation is particularly relevant in practice at high frequencies, where the 
density of resonant frequencies is higher. 

The material is organized as follows: the problem of interest is presented in Section [2j and coupling 
procedures are detailed in Section [3j where the main mathematical results are stated. The finite- 
dimensional approximation of the formulations is addressed in Section [4j Finally, numerical results are 
presented in Section [5j and some conclusions are drawn in Section [6j In Appendix |Aj some details are 
given on the application of the Prandtl-Glauert transformation to the convected Helmholtz equation. 
The proofs of the mathematical results are presented in Appendix [Bj 

2 Aeroacoustic problem 

This section describes the problem of acoustic scattering by a solid object in a non-uniform convective 
flow, together with the underlying physical assumptions. 

2.1 Notation and preliminaries 

Figure [T] describes the geometric setup. The interior domain, corresponding to the area near the 
scattering object where the convective flow is non-uniform, is denoted by Sl~. In the exterior domain, 
the convective flow is assumed to be uniform. The complete medium of propagation, denoted 
by SI C M 3 , is such that SI := f2+ U fT U r ro = M 3 \{solid object}, where r ro := <9fi+ n dQ~ is the 
boundary between the interior and exterior domains. The surface is assumed to be Lipschitz. 
Such an assumption is sufficiently large to include for instance polyhedric surfaces resulting from the 
use of a finite element mesh in Q~ . The surface of the solid scattering object, 9Sl _ \r oo , is denoted by 
T and is assumed to be Lipschitz. 

The speed of sound when the medium of propagation is at rest is denoted by c, the wave number 
by k, the density by p, and the acoustic velocity and pressure, respectively, by v and p. The rescaled 
velocity is defined as M := c _1 v, where M := \M\ is the Mach number. The subscript oo refers to 
uniform flow quantities related to Q + , whereas the subscript refers to point-dependent flow quantities 
related to 0", that is, p\ n - = po(x), p\ n + = p^, k\ n - = k (x), k\ n + = koo, cp- = co(x), C| n + = c^, 
M\q- = Mq(x), M\q+ = Moo- The convective flow is continuous through and tangential on T. 
Hence p, k and M are continuous through T^, and M ■ n = on T. 

The source term g is time-harmonic with pulsation u and is assumed to be located at J7 + for 
simplicity. This source term is an acoustic monopole located in x s € f2 + of amplitude A s , so that 
g := A s 5 Xs cos(wt), where 5 denotes the Dirac mass distribution. The physical quantities are asso- 
ciated with complex quantities with the following convention on, for instance, the acoustic pressure: 
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Figure 1: Geometric setup for the coupled problem 



p <-> Re (pexp (— iut)). In what follows, we always refer to the complex quantity. Furthermore, the 
Hermitian product of two vectors U, W G C 3 is denoted by U ■ W := Yli=i U%Wi, where 7 denotes 
the complex conjugate, and the associated Euclidian norm in C 3 is denoted by || • ||. 



2.2 The convected Helmholtz equation 

In the interior domain Q~ , the convective flow is supposed to be stationary, inviscid, isentropic, 
irrotational and subsonic. The acoustic effects are considered to be a first-order perturbation of this 
flow. With these assumptions, the acoustic velocity v is irrotational, so that there exists an acoustic 
potential <p such that v = V(/j. 

Following [23, Equation (F27)] and [20J, and making use of the irrotationality of the convective 
flow, the linearization of the Euler equations leads to 

p (k 2 (f + ikM ■ V(p) + V • {p(V<p - (M ■ Vip) M + ikipM)) = g in £1, (1) 

where <p is the unknown acoustic potential, and p, k, M, and g = A s 5 Xa are known. Equation 
([!]) is the convected Helmholtz equation. Under the assumption that the acoustic perturbations 
are perfectly reflected by the solid object, the acoustic potential verifies an homogeneous Neumann 
boundary condition on F: 

• n = on T. (2) 



Problem is completed by a Sommerfeld-like boundary condition at infinity (see (12c) below) 



that selects outgoing waves to ensure uniqueness of the solution. In the exterior domain Q + where 
the flow quantities are uniform, equation simplifies into 



A(p + ki (p + 2ik 00 M 00 -V(p-M 00 -V(M 00 -V(p)=g in ft+ (3) 

If there were no scattering object and if the convective flow were uniform in M 3 (and thus equal to 
the flow at infinity), the source term g would create an acoustic potential denoted by ip- mc in M 3 . This 
potential, which solves ^ in M 3 , has an analytical expression, and tp- mc and n ■ V<^i nc are continuous 
across T^. The acoustic potential scattered by the solid object is defined as (p sc := — tp\ nc in 
Eliminating the known acoustic potential ip- mc created by the source yields 

Ai Psc + hl (p sc + 2ik 00 M 00 -Vtp sc - Mco-V {Mco-Vcpsc) = inQ + . (4) 



3 



2.3 The Prandtl— Glauert transformation 



The Prandtl-Glauert transformation was introduced by Glauert in 1928 [19] to study the compressible 
effects of the air on the lift of an airfoil and was applied to subsonic aeroacoustic problems by Amiet 
and Sears in 1970 pQ. Herein, the Prandtl-Glauert transformation is applied in the complete medium 
of propagation and is based on the reduced velocity Moo- This transformation consists in changing 
the space and time variables as 

x' = 7oo (M M • x) Moo + (x- (Moo • x)Moo) x £ tt, 

7 2 ' ' (5) 

t' = t - —Moo x t£R, 

Coo 

where 7oo '■= -, 1 and := M 00 1 Moo with Mqo := |Moo|- The spatial transformation cor- 

V 1 ~ML 

responds to a dilatation along Mqo of magnitude 7oo> the component orthogonal to iVfoo being un- 
changed. In what follows, we suppose that Mqo < 1, so that the Prandtl-Glauert transformation is a 
C°°-diffeomorphism from f2 x M to Vt' x M, where £1' denotes the transformed medium of propagation. 

2.4 The transformed problem 

Let / be such that y>(x) = f(x')exp(—ik 00 'y 00 (M 00 -x')), x' £ Cl'; f- mc and / sc are defined from 
</?inc and ip sc in the same fashion, so that f inc is analytically known, and defined in M 3 . Let s(x') := 
p" 1 exp (ifcooToo (Moo • x')) g(x'), x' G 

In what follows, the transformed geometry, unknowns and operators are considered unless specified 
otherwise. For brevity, primes are omitted. The convected Helmholtz equation ([!]) after the Prandtl- 
Glauert transformation becomes 

rk 2 (3f + irkV • V/ + V • (irkfV + rEVf) = q in O, (6) 

while the boundary condition ^ becomes 

(irkfV + rEVf) ■ n = on T, (7) 

where r := := (1 + qP) 2 - q 2 Ml, V := (1 + qP)MM - q^M^, q : = 7 ^%S P := M ■ 
and E := MOM with M := I + CooMooM^, O := I — MM T , and Coo := The derivation of 

equations ([6])-([7]) is sketched in Appendix [Aj An important observation is that in Q + , (3 = 7^, V = 
and S = /, so that ^ becomes 

A/ + = ? in^+, (8) 
where fcoo := 7oo^oo- Moreover, since supp(?) C /; nc satisfies 

A/inc + klofinc = ? inO+, A/ inc + ^/ inc = inR 3 \f}+. (9) 
Eliminating f mc in Q yields, 

A/bc + fc^/sc = inft+. (10) 

This is the classical Helmholtz equation with modified wave number feoo • Another important property 
is that the matrix E is symmetric positive definite in Cl~ as well. Indeed, there holds in Q~ , 

C7 • SJ7 > (1 - M 2 ) ||£/|| 2 for all U € C 3 , (11) 

where Mo is uniformly bounded away from 1 since the convective flow is assumed to be subsonic. 
Moreover, for all U, W G C 3 , there holds U ■ EW < ^§-\\U\\ [| W\\ in fi". 
In summary, the boundary value problem we consider is 

rh 2 (3f + irkV • V/ + V • (irkfV + rEVf) = ? in 0, (12a) 

(zrfc/V + rSV/) • n = on T, (12b) 

lim r f 8(/ ~ /inc) - ikotf - f mc )) = 0, (12c) 
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where / is searched in H} c (0) := {u G H 1 (K),yK C compact}. Equation (12a) is the transformed 
convected Helmholtz equation, ( |12b[ ) the transformed boundary condition, and the condition at in- 
finity (12c) the classical Sommerfeld radiation condition that guarantees existence and uniqueness 
for Helmholtz exterior problems |31l Theorem 9.10]. In the general case, the Sommerfeld radiation 
condition is written for the scattered potential, since some incident acoustic potentials, e.g., plane 
waves, do not verify it. 



3 Coupling procedure 

Problem ( |12[ ) is separated into an interior problem posed in Q~ and an exterior problem posed in 
£l + in view of using different numerical methods in each subdomain. Specifically, the problem in the 
interior domain is solved by means of finite elements, whereas the problem in the exterior domain is 
solved by means of boundary elements. The main purpose of this section is to derive two coupling 
procedures between the interior and exterior problems. 



3.1 The transmission problem 

The one-sided Dirichlet traces on of a smooth function u in r2 + UO~ are defined as 7o = it =t: 
and the one-sided Neumann traces as J^u^ := (Vm^) Ir^ ■ n, where := u\q± and where n is the 
unit normal vector to conventionally pointing towards Q + (see Figure [l]). The one-sided normal 
traces on of a smooth vector field <x in Sl + U are defined as 7^cr =t := cr^lr^ • n, where := 
ct\q±. These trace operators are extended to bounded linear operators 7^ : H 1 ^^) — > /^(r^), 
jf : H(A, 0±) -> ff-l^oo) and 7 ± : H(div,n ± ) fl-sfToo), where H 1 ^), iT^r^), and 
H~HToo) are the usual Sobolev spaces on 0± and T^, H(div,Q ± ) = {a £ L 2 (0 ± ), V-cr € L 2 (0 ± )}, 



with L 2 (J7 =t ) the Lebesgue space of square integrable functions on fl^, and H(A,^ 



i7 1 (0 ± ), Au^ 1 S L 2 (r2 ± )} (see [371 Lemma 20.2]). It is actually sufficient to consider functional spaces 
on compact subsets of Q + to define exterior traces on Too. Let X denote the surface V or Too. The 
L 2 (X)-inner product (•, •) : L 2 (X) x L 2 (X) ^ C is defined as 



x 



X(y)fi(y)ds(y). 



(13) 



This inner product can be extended to a duality pairing on H 2 (X) x H2 (X). 

We consider the following transmission problem where the one-sided normal trace 7~ r on T from 
is used to formulate the boundary condition ( |12b ): 



lim r 



d(f + ~ g 
Or 



irkfV + rSV/)" 


= 


in Q, , 


(14a) 


Af sc + fsc 


= 


in n + , 


(14b) 


irkfV + rEV/y 


= 


on r, 


(14c) 




= 


on Too, 


(14d) 


7i7 + -7rr 


= 


on Too, 


(14e) 




= 0. 




(14f) 



Proposition 3.1. Problem (12) is equivalent to Problem (14). 



Proof. If / solves (12), it is clear that / solves (14). Conversely, let / be defined in Q + U $7 such 



that / verifies (14). From |31| Lemma 4.19], since / solves (14a) and (14b), then / solves (12a) if and 
only if the jumps on of its traces and of the normal component of irkfV + rSV/ vanish. The 
first condition is just (14d), while the second condition is 7 1 h / + — 7^ (irkfV + rEV/) - = 0, which, 
by continuity of the convective flow across T^, is just (14e). Finally, (12b) and (12c) are simply (14c 



and (14f) 



Theorem 3.2. Problem (14) is well-posed. 



Proof. See Proof B.6 



□ 



□ 
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3.2 Basic ingredients of the coupling procedure 

The coupling procedure hinges on a weak formulation in the interior domain and a Dirichlet-to- 
Neumann map (DtN) associated with the classical Helmholtz equation ( 14b ) in the exterior domain 

n+. 



3.2.1 Derivation of Dirichlet-to-Neumann maps 



For u G i? 1 U ^ _ ), the jump and average of its Dirichlet traces across are defined respectively 

as [70^]^ := 7o~ u+ ~~ 7o~ u_ an d {70^^ := | (7o~ n+ + 1o u ~)- ^or u e ^ + U the jump 

and average of its Neumann traces across are defined respectively as [7iu] r := 7i~"U + — 7f u~ 
and {7iu} Foo := \ (7^+ + ^u~). 

In what follows, Helmholtz equations, as well as corresponding boundary integral operators, are 
written for the transformed wave number k^. A function u defined over R 3 is said to be a piecewise 
Helmholtz solution if u\q+ and u\ 
fi+ and M 3 \0 H 



solve the classical Helmholtz equation (14b) respectively in 



A radiating piecewise Helmholtz solution is a piecewise Helmholtz solution that 

the single-layer potential 



satisfies the Sommerfeld radiation condition (14f). For all A G C°(T 



is defined as 5(A) (x) := J r E(y — x)\(y)ds(y) , x G M 3 \r oo , where E(x) 



cxp(ifc 00 |a3|) 
4tt\x\ 



is the 



satisfying the 



fundamental solution of the classical Helmholtz equation (14b) with wave number koo 
Sommerfeld radiation condition (14f). For all fi G C°(r 
V{v)(x) := f Toc V y E{y - x)n(y)ds(y), x G 

be extended to bounded linear operators S : H'^T^) -> i^ c (M 3 ) and V : ^(r^) i7 1 1 oc (]R 3 \r 
Moreover, both map onto radiating piecewise Helmholtz solutions. Recalling |25| Theorem 3.1.1], a 
radiating piecewise Helmholtz solution u can be represented from its Dirichlet and Neumann jumps 
across in the form 



00;, the double-layer potential is defined as 
3 \r oo . From [35] Theorem 3.1.16], these operators can 

H{ oc (R 3 ) and V : H^T^) -> H}, 



00 j 



u = -S([ 7l «] r J + D([ 7 o«]r J in ft+ U (R 3 \ft+). 

The single-layer and double- layer potentials satisfy the following jump relations across T c 
rem 3.1.2]: 



[70 (SX)] roc = 0, [71 (5A)] roo = -A, VA G H-^(Too), 

[70 (Pfi)] rao = m, [71 CD/i)] r =0, v M g m (roo), 



(15) 
Theo- 



(16) 



with the operators 



(17) 



S-.H-* (roo) -►#3 (Too), 5A:= 7o (5A), 

D : fflCToo) -»• Fl^oo), £>// := {70 (2?/x)} r<>o , 

£> : H-zpoo) -> F-s^oo), L>A := {71 (SA)} roo , 

NrfTsCrooJ-^frs^oo), iV/i := -71 (D M ) , 

being respectively the single-layer, double-layer, transpose (or dual) of the double-layer, and hyper- 
singular boundary integral operators. The Dirichlet and Neumann traces are well-defined, and the 
mapping properties can be found in |31l Theorem 7.1]. The following trace identities are directly 



derived from (16): 



70S = S, 



76^ 



1 



D±-I, 
2 ' 



7^ 
7l D 



^ 2 ' 
-N, 



(18) 



Moreover, if u is a radiating piecewise Helmholtz solution, taking the interior traces of (j 1 5 p and using 

D 



( 18 ) leads to 



N 



S 

\i + b 



7i u 



(19) 



Let now / solve (14) and let v be the function defined by v\q+ := f sc and i>|r3\o+ := — / inc - The 



function v is a radiating piecewise Helmholtz solution (on f2 + this follows from (14b) and (14f), and 
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on M 3 \f2 + from ([9])). Since f inc is continuous across Too, 

hav] Too = lo fsc + 7(7 /mc = 7o" /sc + 7o fL = 7 + / + - 



(20) 



Likewise, [7iu]p = it f + ■ in what follows, we drop ± superscripts for the Dirichlet and Neumann 



traces of / and f- mc since the traces are continuous. Then, (19) applied to u = v yields 



\l-D 
N 



S 

J + D 



70 / 

71 / 



7o /i 
71 fi 



inc 
inc 



(21) 



Various identities relating 70/ and 71/ can be derived from (21 ), and these identities can be used to de- 



fine DtN maps. For example, using the first line of ( 21 ), 71/ = DtNo(jof) := (70/inc + (D — 70/) , 
where the question of the inversibility of S has to be addressed. Two other examples are detailed in 
Sections 13.31 and 13.41 below. 



Remark 3.3. The block operator defined on the left-hand side of (21) is not injective. 



3.2.2 Weak formulation in the interior domain Q 



Let <3? := /|q- where / solves (14). Multiplying ( |14a ) by a test function $* G i7 1 (J7 ) and using a 
Green formula yields 



0. 



(22) 



with the sesquilinear form 



r~V<I> • V$ 



n- Ja- 



~/ nr (irk$V + rEV$) ,7 r $ * 



and the boundary term Q-p :- 

Dirichlet trace on T from Q~ 
domain fT is: Find <S> £ H 1 (fT) such thatW* G H 1 

V($,**)-(7r*>7o~**) r 



$*V$) (23) 

where 7^~ r denotes the one-sided 
Owing to (14c), Qr = 0. Hence, a weak formulation in the interior 



0. 



(24) 



Using the transmission conditions ( 14d)-( 14ej ), 7 $ = 70/ and j x $ = 71/, so that the coupling 
with the exterior problem can be written as 7^$ = DtN^yQ). This yields the following coupled 
formulation: Find <3? G H 1 (fi - ) such that V$* G i? 1 



V(*,$*) - (£>tiV( 7 o $),7 ~$') roo = 0. 



(25) 



3.3 Unstable coupled formulation 



To carry out the coupling, a first classical DtN map is considered. Since this DtN map is not well- 
defined at some frequencies, the resulting coupled formulation is not well-posed at these frequencies, 



and is therefore called unstable. From (21), recalling 7 = 70/ and 7 1 = 71/, there holds 



\I-D 
N 



S 

J + D 



7o_$ 
7i * 



7o/inc 

71 /inc 



(26) 



Using the first line of (26), 7l $ = S 1 ((-D — (70 3>) + 7o/inc)- At this point, the inverse of S 
is written formally. Conditions of inversibility are discussed below. From the second line of (26), 



7-f <& = — iV(7 <£) + — Dj (7 1 $) + 7i/in C - Injecting into the right-hand side of this relation the 

expression of 7j~<J? derived above yields the DtN affine map: -DtiV unsta b ; i75(r oo ) — > ff~5(r oo ) such 
that 



7 *I> - £>/ Anns, ah (77*) :- -A">,7*) + ( \l~ D 1 ,S- 



1 



D - - J (70 + 70/inc + 71 /inc (27) 
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The operator inversion requires to introduce the auxiliary field A G H 2 (Too) such that 

1 



D 



l\( 7 -^)-SX = - l0 fn 



yielding 



^AWstab(7o"<*>) = -N(%<5>) +[-I-D\(X) + 7l /i 



(28) 



(29) 



Injecting DtN unst3 ^(j <&) from (29) into the formulation ( [25] ) yields, using (28), the following coupled 
variational formulation: Find (<I>, A) G H such that, V (<&*, A') G "H, 



(30a) 
(30b) 



V(*,**) + (iV( 7 o-$),7o~^) roo + ( [D- 2 I) (A), 70 -** J = (71 Zinc, 7 "^) roo 

\ \ / / Too 

- 2 7 ) (7o ^ } )r " (Ai ' 5(A))r - = ~ ^' 70/inc )r oo • 
with product space H := tf 1 (tt-)xH-2 (Too) and inner product (($, A) , ($*, A*))~, : 



ifi(n-) 



+ 



(A, A*)^_i . The formulation (30) is called unstable since it admits infinitely many solutions at 
some frequencies of the source, leading to incorrect numerical results. 

Remark 3.4. The -DiiV unsta b affine map was proposed by Costabel to obtain a symmetric coupling in 
the case of self-adjoint operators }10\j . The DtN nnstai h map can be well-defined for certain operators: 
for instance, for transmission problems for the Laplace equation, this map leads to a well-defined 
symmetric system. In the system (30), the only non-symmetric contribution results from the vector 
V in the sesquilinear form V. The system becomes symmetric when the flow is uniform everywhere. 
However, since D and D are dual but not adjoint operators, there is no Hermitian symmetry. 



Proposition 3.5. Iff solves (14), then (/ ,71/) solves (30). Conversely, if(<&,\) solves ( |30[ ), then 
K($,\) solves Q, where K : U -> H\ oc {^AY^) is such that Tl($,X)\ u - := $ and K{$,\)\ n + := 

(-5(A)+2?(7o*) + / in c)|n + . 



Proof. See Proof B.7 



□ 



The main difficulty with the unstable coupled formulation (30) stems from the fact that ker(5 l ) 
depends on whether — k^ belongs to the set A of Dirichlet eigenvalues for the Laplacian on the bounded 
domain M 3 \r2 + . Specifically, ker(S') = {0} if — ^ A, while ker(S') contains nontrivial elements if 

-hie A. 



Proposition 3.6. If f solves ( |14[ ), then for all A* G ker(S), (/ _ ,7i/ + A*) solves ( |30| ). 
Proof. This is a direct consequence of ker(5) = ker(L> — \l). 



□ 



Theorem 3.7. If —k^ £ A, then problem (30) is well-posed. If —k^ G A, then (30) admits infinitely 
many solutions of the form (/ _ ,7i/ + A*), where f is the solution to (14) and A* is any element in 
ker(S). 



Proof. See Proof B.8 



□ 



Remark 3.8. Let —k^ G A. Owing to Proposition 3.5, for any couple ($,A) solving (30), TZ(<&,\) 
solves (14). However, even if our goal is to solve (14), we will see in Section^ that the numerical 
procedure to approximate (30) fails to the point that TZ(&,\) is dominated by numerical errors. 
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3.4 Stable coupled formulation 

The idea of considering a linear combination of S and ^1 + D to derive well-posed boundary integral 



equations was independently proposed in 1965 by Brakhage and Werner [5 J , Leis |27J and Panich 
This is the so-called Combined Field Integral Equation (CFIE). However, S and D map H~z (Too) into 
different spaces (H 2 (Too) and H~z (Too) respectively). This inconsistency in the functional setting can 
be solved by considering a regularizing compact operator from H - ?^^) into H^iT^), as observed 
by Buffa and Hiptmair [7J. We briefly recall the approach of [7] and apply it to the present setting. 
Let Vr^ denote the surfacic gradient on Tqq . Consider the following Hermiticin sesQuilinear form: For 
all ^qeH 1 ^), 

froo (p, q) ■■= (V roo P, V roo g) roo + (p, q) Fao , (31) 

and the regularizing operator M : H~ 1 (F OQ ) — > iif 1 (r oo ) defined through the following implicit relation: 
For sdlpei^CToo), 

<5r 00 (Mp ) g) = (p,g)r 00 , Vgeff 1 ^). (32) 

It is readily seen that M = (— Ar^ + 1) -1 , where Ap^, is the Laplace-Beltrami operator on Too. 

Many choices of DtN maps based on CFIE strategies with the regularizing operator M lead to well- 
posed systems whatever the value of k^. The present choice hinges on the inversion of the operator 
S + ir]M (^I + D^j mapping iJ~2(Foo) into #2 (Too) since, from [7J Lemma 4.1], this operator is 
bijective as long as the coupling parameter 77 is such that Re(r/) ^ 0. To do so, the first line of (26) 
and the application of M to the second line of (26) are used to obtain 



|l — D) + irjMN S + ir]M[\l + b 



N 



ll + D 



7o_$ 
7i $ 



7o/inc + ir}Mjif h 

7l/inc 



(33) 



Then, using both equations in (33) in the same fashion as in Section 3.3 leads to DtN sta b '■ (Too) 
-£^ 5 (Foo) such that 



7f* =DtN stah (%$) := -N(%$) +[ 2 I-D 



S + ir)M [^I + D 



1 



I - D + irjMN 



(l0 $) + 7o/inc + ir} M llfmc + 7l/i 



(34) 



The operator inversion requires to introduce the auxiliary field X £ H 2 (r ra ) such that 



•S'+ ///.U ( ^ + £ 



(A) + 



1 



I - D + irjMN 



( 7o $) = 7o/inc + «W(7l/n 



so that 



^JVstab(7o $ ) 



-iV( 7o -$) +(-!-£») (A) +71/1. 



(35) 



(36) 



The evaluation of M involving an operator inversion as well, it requires to introduce another auxiliary 
field p G H l (Too) such that, for all q G ff 1 (r oo ), 



^r oo 0,9) = (iV(7 -$), 



+ 



1 



/ + £> (A) ,9 



(71 Zinc, ?) 



(37) 



so that equation (35) can be rewritten 



S(X) + Q J - D) ( 7o -S) + iryp = ^/fa 



(38) 
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Injecting -DtA^ sta b(7 3>) from (36) into the formulation (25) yields, using ( [38] ) and (37), the following 
stable coupled variational formulation: Find ($, X,p) £l such that V (<£*, A*, p) G H, 



V(<&, **) + (iV( 7o -$), 7o"^)r + D - o J ( A )^o"^ 



(71 Zinc, 7 **) r 



A*, ( ^ - \I ) (To*) ) " (A*, S(A)) rco - zr? (A*,p) roo = - (A*, 70/inc)^ , 



(iVfro + D + -J (A),p 



(39a) 
(39b) 
(39c) 



with product space HI := H 1 (0 )xi? 2 (r^) x H 1 (Too) and inner product (($, A,p) , ($*, A*,p*)) I 
(*>*V<n-) + ( A ' A Vicr-) + ("%'(r»)' 



Proposition 3.9. /// solves (14), i/ien (/ , 7 i/, 0) solves (39). Conversely, if(&,\,p) solves (39), 



i/ien 72.($, A) solves (14) and p = 0, where 1Z is defined in Proposition 3.5 



Proof. See Proof B.3 



Theorem 3.10. Problem (39) is well-posed at all frequencies. 



Proof. See Proof B.5 



□ 



□ 



Remark 3.11. Hiptmair and Meury \22§ derived abstract trace transformation operators and general- 
ized Calderon projectors for the Helmholtz transmission problem. By construction, integral operators 
written using these projectors enjoy the uniqueness property whatever the value of k^. The map 
DtNgta^ corresponds to a particular choice of the trace transformation operator in the general setting 



of J1H Section 9]. Hence, Theorem 3.10 is an extension of plj, Theorem 7.3] to the non-zero flow 
case, for a particular choice of the trace transformation operator. 



4 Finite-dimensional approximation 



The coupled formulations (30) and (|39h are approximated by finite element and boundary element 



methods. The underlying results are well-known from both theories and can be directly applied to 
the present setting. 



4.1 Discrete finite element spaces 

Let Ai be a shape-regular tetrahedral mesh of 0~ . The mesh Too of r ro is composed of the boundary 
faces of M.. Let Km > denote the mesh size, Vj^ the space of continuous piecewise affine polynomials 
on Ai, Sji^ the space of piecewise constant polynomials on Too, and S]^ the space of continuous 
piecewise affine polynomials on Too- Let Hm := VjU x ^tw; anc ^ ^-M '■= V]U x $m x &m- The 



discretization of (30) reads: Find (&m>^m) £ Hm such that, V (^m^^m) e "H-Mi 

a un S tab ^ Mt X M ) , (& M ,\* M )) = 6 UDStab (& M ,\* M ) , 



(40) 



with a unstab and 5 unstab readily deduced from ( |30[ ), while the discretization of (39) reads: Find 
(®mAm,Pm) G such that, V[®mAm,Pm) g U M, 



,stab 



{($m,*m,pm), (&m,*m,Pm)) 



jstab I ^.t 



(41) 



with a stab and b stab readily deduced from ([39]). Since tLm C tL and HLvi C H, both approximations 
are conforming. 
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In what follows, A < B denotes the inequality A < cB with positive constant c independent of the 
mesh size and of the discrete and exact solutions. The following classical approximation properties 
are available (see [H HH [35] ) : 



inf 11$ - $M\\mfn-\£ h M II* 



inf || A — Ax | 



< h 



M 



Hi (Too) 



inf ||p-PA4||/ji( roo ) ^ ^ H^r.,) • 
Hence, the following approximation properties hold: V(<3?,A) G H 2 (Q~) x /^(r^), 



(42) 



^ inf m,\)-($MAM)\\ H Z h M(\\nHHn-) + \\M\ H h r , 



(43) 



and V(*,A,p) G H 2 (n-) x ^(r^) x iF^), 

, inf , ra ll(^A,p)-($x ! A M ,px)|| H </i M (||$|| i72(n _ ) + ||A|| i ^ , + Wfla ( r 00 ))- ( 44 ) 

Remark 4.1. Taking a polynomial approximation with one order less for H~? (^oo) than for H 1 (fi _ ) 
and H l (r m ) enables all the approximations to be at the same order in hj^\. 



4.2 Discretization of the coupled formulations 

Let (Oi)i<i< p and (*pi)±<i< q denote finite element bases for Vj^ and Sj^ respectively. The decompo- 
sitions of <&m G V]^ and Xm G Sj^ on these bases are written in the form &m = Yli=i ^Mi@i an d 
= Yh=i ^Mi^i- Let 



unstab _ / (®Mi) l<i<p 
' ' " (AMi) l<i<<j 



B 



unstab _ I (7l/inc,7o ^VToo 1<*<P 
I - (l/H, 70/inc) roo i<j<g 



(45) 



t unstab 











^-3 J ) (7o"%)) roo 


-(^,5(^-)) roo 



(46) 



where in J 4 unstab the index i refers to the lines and the index j to the columns. The linear system 



resulting from (40) is 



^junstab^unstab ^unst-ati (47) 

Let (£i)i<i< r denote a finite element basis for Sj^. The decomposition of pm £ $m on this basis 



is written in the form p_M = Yli=iPMi^i- Let 



stab 



{^Mi) l<i<p 
(AMi) l<i<q | ■ 
(PMi) l<i<r 



B 



stab 



( (71 Zinc, 7 ^)r,» l<i<p 

~ (^i,70/inc) roo \<i<q I > 
y (7l/ino l<j<r 



(48) 







((D-ij) Wi).7oft) rm 





^stab 


(i>i,(D-±I) (7 -%)) r ^ 










((S-J/)(W,6) r _ 





(49) 



/ 



with the same convention on the indices i and j of yl stab . The linear system resulting from (41) is 



4 stab stab nstab 



(50) 
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4.3 Inf-sup stability of the discretized formulations 



From the Fredholm setting and the approximation properties (43) and (44), the following discrete 
inf-sup conditions can be derived following |24} Theorem 14]. 

Proposition 4.2. If —k^ ^ A and h_M is small enough, there holds, for all ($m-i^m) £ T~Lm> 



sup 

(0,0)^(4 t M ,A t M )6«A1 



II ($M, Am) || 



At all frequencies and if h_M is small enough, there holds, for all (^m^M^Pm) S 1^, 

|a stab (($mAMiPm),(®m,X M ,Pm))\ 



sup 

M 



^ ||(^M,Ax,PMy 



(51) 



(52) 



4.4 Convergence 

From the inf-sup stability of the discrete problems, the following error estimates easily follow from 
[241 Theorem 13]. 

Proposition 4.3. If —k^ A and Hm is small enough, the discrete problem (40) has a unique 
solution (<&Mi^m) S T~Lmj an d the following optimal error estimate holds: 

\\{®,\)-{®mAm)\\h< , inf \\{*,*)-{*m,>?m)\\h, (53) 



where (<!>, A) is the unique solution of (30). At all frequencies and if Km is small enough, the discrete 
problem (41) has a unique solution ($^, ^MiPm) £ HU/j, and the following optimal error estimate 
holds: 

\\($,X,p)-($M,*M,Pm)\\m< , inf ||(*,A,p)- ($5wA>PA"l) Ik (54) 

(*5w. A 5w.P$w)eliM 

where (<£,A,p) is £/ie unique solution of (39). 

Remark 4.4. T/ie constant in ( |53[ ) depends on k^, and its value explodes as —k^ tends to an 
element of A. The constant in (|54l) depends on /uqo ^s well, but remains bounded on any bounded set 
of frequencies. 



5 Numerical results 



The purpose of this section is the comparison between the unstable formulation (30) and the stable 



formulation (39) with the coupling parameter 77 = 1. Both formulations have been implemented in the 
EADS in- house boundary element software called ACTIPOLE |12|ll3j. This software can treat general 
three-dimensional geometries. The iterative solver is a block generalized conjugate residual method 
[32, 26] based on a generalized minimal residual method [34], suitable for non-symmetric systems. A 
sparse approximate inverse preconditioner JBJ [9] is used. 

Consider an ellipsoid with major axis directed along the z-axis. This object is included inside a 
larger ball. The external border of the ball after discretization is the surface Too. A potential flow is 
computed around the ellipsoid and inside the ball, such that the flow is uniform outside the ball, of 
Mach number 0.3 and directed along the z-axis. An acoustic monopole source lies upstream of the 
ball, on the z-axis as well. Four different meshes are considered, see Table [T] Mesh 1 is such that the 
mean edge is ten times smaller than the wavelength at rest for a frequency of 1500 Hz. From Table 
[TJ for fine meshes, the number of basis functions used to discretize the unknown p for the variational 



formulation (41) takes a smaller part in the total number of basis functions than for coarse meshes. 



Therefore, the relative complexity added to (30) by the third equation of (39) decreases with the total 
number of unknowns, which is an interesting property when it comes to industrial test cases. Figure 
[2] displays Mesh 1 and the rescaled velocity Mq of the potential flow. 

In what follows, a frequency f is called resonant if — = — 4?r 2 f G A, where A is the set of 

Dirichlet eigenvalues for the Laplacian on M 3 \fi + . The set A depends on the geometric shape of the 
coupling surface T^, which slightly changes after each discretization. 
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Mesh 1 


Mesh 2 


Mesh 3 


Mesh 4 


number of volumic dofs <3? 


1796 


687 


194 


79 


number of surfacic Po dofs A 


808 


510 


270 


148 


number of surfacic Pi dofs p 


406 


257 


137 


76 


proportion of dofs p in the total number of dofs 


11.9% 


15.0% 


18.6% 


20.0% 


smallest edge (mm) 


7.09 


8.78 


15.71 


19.18 


mean edge (mm) 


22.64 


32.20 


49.78 


66.46 


largest edge (mm) 


56.87 


70.62 


103.59 


112.71 




5.1 Comparison of pressure fields 



As seen in Theorem 3.7 the unstable formulation (30) is not well-posed at resonant frequencies. First, 



a prospective study to identify a resonant frequency for each of the four meshes is carried out by 
monitoring the condition number of the matrices produced by the discretized version of the unstable 



formulation (30). A resonant frequency for Mesh 1, Mesh 2, Mesh 3, and Mesh 4 is identified around 



1509.849 Hz, 1513.431 Hz, 1521.015 Hz, and 1535.704 Hz respectively. 

The convergence of the iterative solver is monitored by requiring that the Euclidian norm of the 
relative residual is smaller than 10~ 6 . Additional tests indicate that the discretized solution to the 
stable formulation does not change much below this value of the relative residual. For Mesh 1, away 
from a resonance, say at 1500 Hz, the scattered pressure fields computed with the unstable and stable 
formulations are very similar. This holds as well for the total pressure fields, see Figure [3| At the 
resonant frequency 1509.849 Hz, the unstable formulation (30) yields pressure maps quite different 



from the ones at 1500 Hz, whereas the stable formulation (39) yields pressure maps very similar to 



the ones at 1500 Hz, see Figure |4j From Proposition 3.5 at a resonant frequency, the solutions to 
( 30 ) that differ from the solution to ( [39] ) should not affect the scattered pressure field in the exterior 
domain. Actually, the distortion of the scattered field with the unstable formulation ( 30 ) is the result 



of the significant magnification of discretization and numerical errors by the ill-conditioning of the 



linear system approximating (30). 
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Total pressure 



Total pressure 





Scattered pressure 
102 




Scattered pressure 
102 




Figure 3: Mesh 1, 1500 Hz. Top : real part of the total pressure; left: unstable formulation (30), right: 



stable formulation (39). Bottom : real part of the scattered pressure; left: unstable formulation (30), 



right: stable formulation (39). At this non-resonant frequency, both formulations yield similar results. 



5.2 Auxiliary variable p 

In Figure [5j the left plot indicates that with Mesh 1, the magnitude of p is around 0.5% of the 
scattered pressure. The right plot shows the behaviour of the magnitude of p with respect to the 
stopping criterion of the iterative solver for the four meshes. The finer the mesh, the smaller the 



auxiliary variable p, which is consistent with the fact that the p-component of the solution to (39) 



vanishes (see Proposition 3.9). 



5.3 Comparison of condition numbers 

Figure [6] presents the condition numbers of the matrices resulting from the formulations ( 30 ) and ( 39 ) 
with respect to the frequency. In the left plot, the curves are centered at the resonant frequencies. The 
finer the mesh, the higher the condition number explodes. The width of the peak at the resonance 
does not appear to depend on the mesh. In the right plot, a larger bandwidth is considered with 
Mesh 2. Owing to the frequency sampling (every 5 Hz), some resonances may be missed and the local 
maxima may not be accurately reached (in particular, from the left plot, the local maximum of 7.2 for 



log(cond(M)) at 1513.431 Hz is very underestimated). The stable formulation (39) produces somewhat 
larger condition numbers for the large majority of the frequencies, but, unlike the unstable formulation 



(|30|), it presents no resonance. Moreover, from the Weyl formula, the number of resonant frequencies 

3 



smaller than f increases as f 2 , making the need for a stable formulation even more important for 
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Total pressure 



Total pressure 





Scattered pressure 




Scattered pressure 
102 




Figure 4: Mesh 1, 1509.849 Hz. Top : real part of the total pressure; left: unstable formulation (30), 



right: stable formulation (39). Bottom : real part of the scattered pressure; left: unstable formulation 
(30), right: stable formulation (39). At this resonant frequency, the two formulations yield different 



results. 



simulations at higher frequencies. 



5.4 Convergence 



To further study the impact of the ill-conditioning of the unstable formulation (30) on the computed 



solution, the preconditioning is not used in what follows. First, the value of the acoustic pressure on 
a network of 10000 points located further than 0.5 m from the center of the sphere (therefore in + ) 



is computed using the stable formulation (39) with Mesh 1 at the resonant frequency 1509.849 Hz. 



This computed acoustic pressure is called the accurate pressure. Next, the acoustic pressure on the 
same network of points is computed for different values of the number of iterations of the solver, using 



the unstable formulation (30) and the stable formulation (39) with Mesh 1 at the resonance 1509.849 
Hz. The relative difference between the computed pressure and the accurate pressure in Euclidian 
norm is called the relative error. Figure [7] presents the relative residual and the relative error with 
respect to the number of iterations. With the unstable formulation (30), the relative residual decreases 



irregularly. In particular, it stays constant during around 200 iterations. The relative error decreases, 
stays constant, rises after 400 iterations, and finally stabilizes at a large value, whereas the relative 
residual keeps converging to zero. As for every ill-conditioned problem, the relative residual cannot 
be used to ascertain convergence towards the correct solution. In particular, the relative residual is 



extremely small, while the error is of order one. With the stable formulation (39), the relative residual 
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Figure 5: Stable formulation (39) at 1500 Hz. Left : real part of the scattered pressure and auxiliary 
variable p with Mesh 1. Right : Magnitude of the auxiliary variable p as a function of the stopping 
criterion of the iterative linear solver with all meshes. 
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Figure 6: Condition number of the matrix for the unstable formulation (30) (solid) and the stable 
formulation (39) (dashed). Left: centered representation around a resonant frequency for the four 
meshes. Right: larger bandwidth with Mesh 2. 



and the relative error decrease regularly, and in the same fashion. 



5.5 Choice of the coupling parameter r] 



In the stable formulation ( 39 ) , the choice of the coupling parameter rj is expected to have a direct 
effect of the condition number of the matrix ^4 stab . In Figure [8j this condition number is plotted 
for Mesh 4 and for various values of rj. For rj = 0, equations (39a)-(39b) are decoupled from (39c), 



and (39a)-(39b) become (30), so that the curve for rj = 0.001 is similar to the curve of the unstable 



formulation for Mesh 4 in Figure |6| The condition number appears to be the smallest for r\ in the 
range 1 to 10, and worsens for lower and higher values of rj. This motivates the choice r\ = \ made in 
the above simulations. 
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Figure 7: Mesh 1 at resonance 1509.849 Hz; relative residual (left) and relative error (right) with 
respect to the number of iterations. 




Figure 8: Condition number of the matrix for the stable formulation (39) centered around the resonant 
frequency at 1535.704 Hz for Mesh 4. In this case, the chosen value rj = 1 leads to the minimal condition 
numbers. 

6 Conclusion 

In this work, we derived two coupled boundary element / finite element methods for the convected 
Helmholtz equation with non-uniform flow in a bounded domain. The first one leads to an unstable 
formulation, while the second one leads to a stable formulation. The unstable formulation involves 
two equations and is well-posed except at some resonant frequencies of the source, while the stable 
one is unconditionally well-posed, but involves three equations. Even if the unstable formulation 
admits infinitely many solutions at resonant frequencies, the pressure field resulting from any of these 
solutions equals the one resulting from the stable formulation. However, our numerical results show 
that at resonant frequencies, the discretization of the unstable formulation is so ill-conditioned that 
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the pressure field is very different from the one produced by the stable formulation. Moreover, the 
stable formulation remains tractable within large industrial problems since the relative complexity 
added by its third equation decreases with the size of the mesh. Its interest is also enhanced by the 
fact that, at higher frequencies, the density of resonant frequencies is more important. 

As long as the uniform flow assumption in the exterior domain is reasonable, more complex flows 
in the interior domain can be considered, as well as more complex boundary conditions at the surface 
of the scattering object. These extensions only require to modify the finite element part of the present 
methodology. An important development for practical simulations is the introduction of modal sources 
in the interior domain, which is the purpose of [2]. 

Another interesting extension of this work is the resolution of parametrized aeroacoustic problems, 
with the frequency of the source as a parameter, using reduced-order models, for instance by means 
of Proper Generalized Decomposition or Reduced Basis methods. Using the unstable formulation 
may involve ill-conditioned numerical resolutions if the frequency range of interest contains resonant 
frequencies, whereas the stable formulation guarantees well-posedness of the procedure. Moreover, the 
complexity of the online stage of the reduced-order model is not increased by the third equation of 
the stable formulation. 
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A Derivation of the transformed convected Helmholtz equation (12a)- 

In this appendix, the geometry, unknowns and operators after the Prandtl-Glauert transformation are 
indicated with primes. To apply this transformation to a PDE in the frequency domain, one has first 
to change the differential operators as 

Vu = MVu, V • U = V • MU, (55) 

for a scalar- valued function u and a vector- valued function U. Here, M = I + CooM^M^ with 
Coo = an d 7oo = 1 M2 > and V refers to derivatives with respect to the transformed variables 

x' . Moreover, it is readily verified that 

MM = M + C 00 PM 00 , MM 00 = loo M 00 , MM ■ = MM^ ■ M = 7oo P, (56) 

where P = M ■ M^. After changing the differential operators, one has to change the unknown 
function as (p(x) = a(x')f(x'), where a(x') := exp (— ifcooToo (-^oo ■ x '))- 

We now show that, following the Prandtl-Glauert transformation, Equation ([!]) becomes 

rk 2 fif + irkV ■ Vf + V • (irkfV + rSV'/) = (57) 

where r = /3 = (1 + qP) 2 - q 2 Ml, V = (1 + qP)NM - q^M^, q = 7 ^%S H = M(I - 

MM T )M, and ? = p^a^g. Dividing Equation by leads to aq = rk 2 (p + irkM • V'cp + V • 
(r (V'cp - (M • VV) M + ikipM)). Applying ((55}, it is inferred 



aq = rk 2 tp + irkM ■ MV'ip + V' • (rMMVtp) - V ■ (r (M • NV'y) MM) + V' • (irkyNM) 

Substituting tp for af and expanding the derivatives with respect to a yields 

aq = ark 2 f + airkM • MV f + arkk^^fiM • MM^) + V' • {arNNV f) 
- V' • (airfc oo7oo /AAAAM 00 ) - V' • (ar (M • MV f) MM) 
+ V' • (mrfcooToo/ (M • MM^MM) + V' • (airkfMM) , 
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since Va = —oik^ 7oe>Afoo- Using (56) and simplifying some terms leads to 

as = ark 2 f + airkM ■ MVf + ark 2 qPf + V' • (arMMV f) - V • (airk^qfM^) 
- V • (ar (M • MVf) MM) + V' • (airk(l + qP)fMM) . 

Expanding again the derivatives with respect to a yields 
C = rk 2 f + irfcM • MVf + rA; 2 gP/ + V' • (rMMV f) - irkoo^ooMMVf ■ 
- V ■ (irkq^fM^) - rkk^lfMl - V' ■ (r (M • MV f) MM) 

+ irk^joo (M ■ MVf) MM • Moo + V • (»rfc(l + qP)fMM) + rA;fc oo7oo (l + qP)fMM ■ M^. 
Using again ( |56| ) as well as the symmetry of M, it is inferred 

<, = rk 2 f + irkMM • V'f + rA: 2 gP/ + V' • (rMMV f) - irkq^M^ ■ Vf 
- V ■ (irfcgWMoo) - rk 2 q 2 M 2 00 f - V ■ (r (MMM T M) V f) 
+ irkqPMM ■ Vf + V • (irk{l + qP)fMM) + rk 2 qP(l + qP)f. 
The terms are now reorganized with respect to the orders of derivation of / to obtain 

c = rk 2 (l +qP- q 2 Ml + qP(l + qP))f 
+ irk ((1 + qP)MM - (rrooMoo) • V/ 
+ V • (frfcf ((1 + qP)MM - (HooMoo)) 
+ V • (rMMV f) - V ■ (r (MMM T M) V f) , 



(58) 



(59) 



yielding (57). 



We now show that, following the Prandtl-Glauert transformation, the boundary condition ^ 
becomes 

(irkfV + rZVf) ■ n' = on V , (60) 

where T' denotes the transformed boundary T. The normals on the initial geometry are denoted by 
n, and the normals on the transformed geometry by n' . It is readily seen that 

n = KooMn' on T, (61) 

where is a normalization factor that is not needed in what follows. Owing to (55) and ( |61[ ), ^ 
becomes MVy-Mn' = 0. Hence, MV(af)-Mn' = 0, leading to (MVf - ik 00 'y 00 fAfM 00 )-Afn' = 0. 
Since the flow is tangential on T, M ■ n = on V. Hence, M ■ Mn' = on T, so that 

k n 



\MVf - (MM ■ V'f) M + ikf ( (1 + qP)M 



-^ loo MM c 



■M' 



n 



(62) 



Using the symmetry of M and ( [56] ) , ( [62] ) leads to ( 60 ) 



B Proofs 



As preliminary results, we first derive Propositions B.l and B.2 Then, recalling that (14) is the 



original transmission problem, (30) the unstable formulation, and (39) the stable formulation, the 



mathematical results are obtained in the following order: 



Proposition B.2 



Proposition B.l 



Link between (14) and (39) (Proposition 3.9) 



Link between (14) and (30) (Proposition 3.5) 



Uniqueness for (39) =^ Well-posedness of p9j) (Theorem 3.10) 



Well-posedness of (14) (Theorem 3.2) 



Conditional uniqueness for ( 30 ) Conditional well-posedness of p0 



(Theorem 3.7) 
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Notice that we do no start with the well-posedness of (14) since, to our knowledge, a result directly 



applicable to it is not available in the literature. Herein, we prove the well-posedness of (39) which is 



equivalent to (14). 



Proposition B.l. Problem (14) has at most one solution in Hh (fi). 



Proof. Let / G H^ oc (Q) solve (14) with q = (so that f inc = 0). From Proposition 3.1, / solves (12) 
Let B be a ball containing Vt~ . Let /* G H(A,B). Using Green's first identity, 



/ (~rk 2 (3f - ikrV ■ Vf - V • (irkfV + (rHV/)) f 
JnnB 



rEVf ■ Vf - rk 2 (3ff - ikrV ■ (Vff - V/*/) - hi, 9B f,lo8Bf\ . „ 



(63) 



where 7 QB and r ) l QB are the Dirichlet and Neumann traces on dB from B. Taking f = f yields 



Tx^b/^b/L = / rHV ^ • V ^ - rk2 Pff - 2krV ■ ( ImV //) 



(64) 



so that Im 9B f, 7 g B fj^ = 0. Using Rellich Lemma (see [3T1 Lemma 9.9]), since / G H^ oc (M. 3 \B) 
solves the classical Helmholtz equation in M. 3 \B and satisfies the Sommerfeld radiation condition, as 



(65) 



well as Im g^/Woae/J^ — ®> ^ s inferred that /|rs\b = 0. Equation Q can be written 



L(f) ■= (rk 2 f3 + V • (irkV)) f + 2irkV • V/ + V • (rSV/) = in 0. 



From |18| Theorem 1.1], since rH is uniformly elliptic with Lipschitz continuous coefficients, and 
rk 2 j3 + V- (irkV) and 2irkV have bounded coefficients, the differential operator L satisfies the strong 
unique continuation property in Q. Hence, /|k3\b = implies that the only H\ oc (Q) solution of (12) 
with ? = 0is_f = 0inf2. The assertion follows from Proposition 3.1 □ 



Proposition B.2. Let (<]?, X,p) solve ( |39| ). There holds p = 0, A = 7 X <£, and 

rk 2 f3$~ + irkV ■ V<I>~ + V • [irk^'V + rHV<J>~) =0 in L 2 (fT), 

7~ r {irk® V + rHV$) = in#~5(r), 

iV( 7o -$) + f £> + J i) (A) = 7 i/ inc in #-3(1^; 



L> - ij) ( 7o -$) - 5(A) = -7o/inc in fI(F x ). 



(66a) 
(66b) 

(66c) 
(66d) 



Proof. Set <r := zr/c^V + rHV^. Owing to (39a) and recalling the definition (23) of V, there holds, 



rk 2 p<5>¥ + i rkV • V$$* - (^(7 ~$) + ^ - -I ] (A) - 7i/inc,7 " $ 



1 



(67) 



Restricting the test function to C£°(f2 ) shows that (66a) holds in L 2 (f2 ). Moreover, since 

^ W-V¥ = - jf _ (V • a) ** + (7-^,70"*%^ + (\r^ 7 >^) r , (68) 
and recalling that 7~<r = 7^$ on Too, owing to the property of the convective flow at Too, 

(aT( 7o -$) + (d - ±l) (A) - 7i/inc,7o~^) + (71^,70 ^) roc + (7^, 7 >*') r = 0- (69) 

\ \ / / Too 
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Then, by the surjectivity of the trace operators 7 from ) onto (Too) and 7 r from H 1 ^ 

i 

onto H? (r) (see 



Theorem 2.6.11]), it is deduced from (69) that (66b) holds, and from (69) and 
( 39b ) respectively that 



A> ; ,7<I>) + ( £> - 1/ J (A) + - a -<j> - - j/u 



J ( 7o -$)-5(A)-^p = -7o/ il 



(70a) 
(70b) 



Owing to (39c) and (70a), it is inferred that 6r ca (p,p t ) = (A — 7 X ^,p*)roo for all p l £ iJ 1 (fi ). From 
the definition Q of M, there holds p = M(A - 7 f$). Let x := A - 7 jf $. Then, from ( |i~9j ), ( ffo| ) can 
be written 



iV 



s 

ll + D 




-irjMx + 70/h 

S + 71 Zinc 



(71) 



so that (—irjMx + 7 o/mc, x + 7i/inc) belongs to the range of the block operator defined in the left-hand 
side of ( |71[ ). Under this condition, from |22[ Theorem 4.1], citing [39], a radiating piecewise Helmholtz 
solution u such that ^qU = —irjMx + jofi nc and 7{~u = x + 7i/i nc can be constructed. Consider v such 
that v\q+ = and u|r3\q+ = u\q- , and w such that w|n+ = and u>|r3\q+ = fine- Since v and w are 
radiating piecewise Helmholtz solutions, u := v — w is also a radiating piecewise Helmholtz solution. 
Since [70 it] = irjMx and [71 it] 



-x. 



D 



(19) implies 



N 



S 

J + D 



irjMx 



-x 




(72) 



Consider now u := S(x) + D(irjMx). From [25, p. 113], the single-layer and double-layer potentials 
are radiating piecewise Helmholtz solutions. In particular, it|]R3\o+ solves the Helmholtz equation in 
R 3 \ft + , therefore 

(7r*T,7o~«~)r =/ {|Vn| 2 -A4|u| 2 }GlR. (73) 



From the trace relation s (pL8[ ) , there holds 7J" vT = 7^ (S(x)+T>(irjMx)) 
using the first line of (|72[), 7p~it~ = —irjMx. Likewise, 7]" u~ = x. From (73), 
However, since (x,Mx) T ^ = \\Mx\\jj 

x = 0, and p = Mx = 0, leading to A 
( |70a[ )-( |70bp using A = 7-f $ and p = 0. 



) 



° and Re(ry) 7^ 0, there holds (x, Mx) T 



Sx+_irj [D - 1) Mx. Then, 
irj(x,Mx) Tao G R. 
■ 0. Therefore 



Finally, (66c)-(66d) are directly obtained from 

□ 



Proof B.3 (Proposition 3.9). 



The first part of the proposition is clear from the results of Section [3] Let (3>,A,p) solve (39). 
Since the single-layer and double-layer potentials are radiating piecewise Helmholtz solutions, TZ(&, A) 
satisfies ( |14b|) and ( |14f|). From Proposition B.2, p = 0, A = 7^ 3>, and (|66j) hold s. (|14a[ ) and ( |14c| ) 
are just (66aj) and (66b). Then, using the definition of 1Z given in Proposition 3.5 and the trace 
identities (jl8), the exterior Dirichlet trace of K(§, A) is JqTZ(<^>, A) = 7^" (-5 (A) + T>(%$) + / inc ) = 
-5(A) + (D+ij)(7^$)+7o/inc- Using ^66dj), there holds 7^($, A) = 7^$, which is 7^^($, A), from 



which we infer that the first transmission condition (14d) holds. The second transmission condition 



(14c) is obtained in the same fashion from the exterior Neumann trace of 1Z(&, A), ( |66c[ ), and using 
the fact that A = 7J" <1>. 



Proposition B.4. Problem (39) has at most one solution 



Proof. Let (<3?,A,p) solve (39) wit h 70/inc = and 71/mc = 0. From Proposition B.2, p = and 
A = 7^<F From Proposition 3.9, 7£(3>,A) solves (14). From the mapping properties of 5 and T> 



(Section 3.2.1), K($, A) e fl£ c (fl + U O 



ft(S,A) G flj c (fi). Then, Proposition 
$ = A)| n - = 0, and A = 7^$ = 0. 



B.l 



Then, from the transmission conditions (14d) and (14e), 
implies that TZ(&, A) = in U. As a result, there holds 

□ 
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Proof B.5 (Theorem [3l0 ). 

Consider the two sesquilinear forms a\ and a 2 on H x H such that 

((<M,P), (<*>', ■= I rHV¥-V^+(iV (7 -$),7o"^)r + W^WhI+ fr.fop') 

Jn- °° 

+ (( 6 °-^^^') r .-( A '( c0 -^)H r j 

a 2 (($,A,p),($*,A*,p*)) :=- / rA; 2 /3¥$* + i / rfcV • ($V$* - $*V$) + ((iV - N°) 7(7$')r 

Jn- Jn- 

+(A',(S-S°)(A)) roo + ( (D - D°) ( A) , ) - (A*, (D - D°) (7o"^)) roo - ivW^)Z 



D + -I)(\)J 



(74) 

where S°, D°, D° and iV° are the boundary integral operators S, D, D and N for koo = 0. Consider 
the linear form b on HI such that 



6($',A*,p*) := (7l/inc,7 + ( A *'70/inc) roo - (7l/inc,P*) r • 



Problem ([39]) can then be written: Find (<&, A, p) G H, such that V ($*, A*,p*) € HI, 

a(($ 5 A,p),(^,AV)) =6($*,A*,p t ), 



(75) 



(76) 



where a := ai + a 2 . We show the well-posedness of ( |76| ) by proving successively that (i) oi and a 2 are 
bounded on HI x EI and b is bounded on HI, (ii) a\ is H-coercive, (hi) the linear map associated with 
a 2 is compact from HI into HI. With the uniqueness of the solution (Proposition B.4), the assertion 
follows by the Fredholm alternative. 



(i) From Section 



2.4 



~V$ • V$< < £tt!M|V$||||V$*||. Then, 
1 



n- 



< 



1 - Ml, 



The other volumic integrals are simply controlled by 
rkV ■ ($V$* - $*V<& 



rk 2 p<$>& 



< 2||rAj|| £ , o( n -)||V|| ioo(n _j3||$|| fr x (:n -)||* Hh 1 ^-)) 

< Hrfc^lUoo^-) ||*||ffi(n-) 11**11^(0-). 



(78) 



From |31l Theorem 6.11], all the involved integral operators are bounded in their natural trace spaces. 
The boundedness constant of an operator A is denoted by Ca, and the continuity constant of the 
interior Dirichlet trace operator is denoted by C^-: \\%^\\ H ^ r ^ — C 7 ~ ll^lli^H^ - )' Moreover, since 

i? 1 (r oo ) C -fT5(r oo ), there exists a constant Cr^ > such that ||p|| i < ||p||#i( roo ). These 

(Too) 

inequalities lead to 
|oi((*,A J jj),(* t ,A t ,p*))| <2 J + 



1 — M 2 



\L°°(n-)\\l + M o\\L°°(n-) + C s o 



+C (1 + C D o + C 6o ) + C 2 C N o || (<&, X,p) || H || (*', A*,j^ 



To 



7o 



|a 2 ((*,A J jj),(* t ,A*,p t ))| <2 



I^ 2 /?IIl-(c-) + 2||rfe|| i cx> (n - ) ||F|| ioc{n _ ) 3 + C roo Q + |r?| + C^) + C s o 



+C S + C 7 - {C Vtx C N + C D o + C D + C 6o + C 3 ) + C 2 ^ (C N o + C N )\ || (*, X,p) || H || (**, A*,p*) 



7o 



fc(^*,A*,p*)| < V2 ((C 7 - +C roo )||7i/inc|l H _i + ll7o/i 



rmc|l ^(r») 



(**,A* 5 p*) 
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(ii) Unlike D and D, the operators D° and D° are real-valued. They are therefore adjoint, so that 

1 



D v --I (A), To"* 



A. ( £>" - ijj ( 7o -$) J e/Ji 



(79) 



From Theorem 2], the operators N° and S° are strongly elliptic in their natural trace spaces. 
Moreover, from Section 



2.4 



for all U G 



U-EU > 1 



M 2 ) || £7" Then, there holds 



Re(oi((A,*,p),(A,* J p)))= / rHV$-V$+(iV (7 -^),7o"*)r + (^S°(X)) Foo + S Voo (p, P ) 
> inf (r (1 -M 2 )) ||V$||ra fn -i + K N o\\%<5>\\ 2 i + -^so||A|| 2 i + \\pfm{r ), 



(80) 



where the coercivity constant of an operator A is denoted by Ka- From the Petree-Tartar Lemma 
|16| Lemma A. 38], it can be shown that there exists a constant Cq- > such that ||3>||ifi(o,-) < 



Cc 



f 2(r 



. Therefore, 



|v$|| L 2 (n - )3 + ho 

Re ( ai ((A,*,p) , (A,$,p))) >^2-min (inf (r (l - M 2 )) H^ 2 ^ 



> min 




ff-*(r„) + W Woo) 



;i-m 2 )) k 



N° 



2CI 



^1- 



K so ,l 



(a-) 



(<M,p) 



(81) 



(hi) Let V be a Hilbert space, let A be an operator from V to V, and let a be a sesquilinear form 
such that, for all u,v £ V, a(u, v) = (Au, v)y. A classical result states that A is compact if and only if, 
for all weakly convergent sequences (u n ), (v n ) 6 V N such that u n — u and v n v, there holds, up to 
subsequences, a(u n ,v n ) ->• a(u,v). Let ($ ln , A in ,pi J^($i, Ai,pi) and (<I> 2n , A 2n ,pi n )^($ 2 , A 2 ,pi) 
be two weakly convergent sequences in H. Since the injection of i7 1 (i7 _ ) into L 2 (S7~) is compact, then, 
up to subsequences, <3?i n — and V<3?j n ^V<]?i, z = 1,2, in L 2 (f2~). Therefore, up to subsequences, 



/ rk 2 (3<5> ln <5>2n + » / • (* in V$ 2n - $ 2 „ V$ in ) 



rfc 2 /3<I>i<]? 2 



(82) 



r/cF • $iVf 2 - $ 2 V<I>i 



Moreover, from [351 Lemma 3.9.8], S-S°, N-N°,D-D° and D-D° are compact operators in their 
natural trace spaces. Hence, up to subsequences, (S — S°) (Ai n )— > (S — S°) (Ai) in (Too), and the 
corresponding results hold for the other boundary integral operators. Since the injection of H 1 (T 00 ) in 
Hz (Too) is compact, then, up to subsequences, Pi n -^Pi in Hz (Too), so that (A 2n ,pi n ) roo — > (A2,Pi) r . 
The last two terms converge as well by continuity of N, 7^ and D, as well as the strong convergence 
of p2 n in Hz (T^) up to subsequences. This concludes the proof. 



Proof B.6 (Theorem 3.2). 



This a direct consequence of Proposition 3.9 and Theorem 3.10 



Proof B.7 (Proposition 3.5). 

The first part of the proposition is clear from the results of Section [3} Let (<&,A) solve (30). In the 
same fashion as in the proof of Proposition B.2, (14a) and (14c) hold, as well as 



\I + D 

-N 





(83) 
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so that (7^"$,7^~$) belongs to the range of the block operator defined in the left-hand side of (83). 
Under this condition, from [22, Theorem 4.1], citing [39], a radiating piecewise Helmholtz solution 
u such that ^qU = and j^u = 7 1 ~<I> can be constructed. Consider the function v defined as 
v\q+ = u\q+ and v|r3\q+ = 0. The function v is still a radiating piecewise Helmholtz solution, and 



its jumps of traces are [70^]^ = 7o ^ an< ^ [ti^Tcx, = 7i ^- From ( fl9j ), (gJ — D) 7 <I> + $ = 
— 7^u~ = 0. Together with the first line of (83), it is deduced that A — 7f G ker(5) = ker(Z) — 



Therefore, (66c) and (66d) hold. Then, since single-layer and double-layer potentials are radiating 



of 7£(<&, A), the transmission conditions (14d) and (14e) are direc tly o btained from (66c 

\l), in the same fashion as in Proof 



piecewise Helmholtz solutions, 7£.(<3?,A) satisfies (14b) and (14f). Finally, taking the exterior traces 

'($, A), the tn 
A-7f $ G ker(S) 



ker(L> 



B.3 



(66d) and 



Proof B.8 (Theorem 3.7) 



Let ($,A) sol ve (p0[ ) with 70/i 
seen in Proof 



B.4 



and 71 



As a consequence, <3? 
ker(S) = {0}, leading to A 

then obtained using the Fredholm alternative by proceeding similarly to Proof 

Let A* G ker(S') = ker(Z) — |/), and / be the solution to (14). From Proposition |3.6| (f~,Jif + A 



0. From Proposition [3^5] Tl (®, A) solves ([14]). As 

that 1l(<f>, A) 



B.l 



A — 7j $ G ker(S). Suppose 



in 11 
fc*, ^ A. Then, 



7?.($, A) G Hy oc (Q). It is then ded uced from Proposition 
rc(*,A)J n - = 

7 1 $ = 0, so that problem (30) has at most one s olutio n; well-posedness is 

Sup pose —hi. 



B.5 



G A. 



solves (30). 
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